Phase oscillators with global sinusoidal coupling evolve by Mobius group action 
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Systems of A'^ identical phase oscillators with global sinusoidal coupling arc known to display 
low-dimensional dynamics. Although this phenomenon was first observed about 20 years ago, its 
underlying cause has remained a puzzle. Here wc expose the structure working behind the scenes of 
these systems, by proving that the governing equations are generated by the action of the Mobius 
group, a three-parameter subgroup of fractional linear transformations that map the unit disc to 
itself. When there are no auxiliary state variables, the group action partitions the A'^-dimensional 
state space into three-dimensional invariant manifolds (the group orbits). The N — 3 constants 
of motion associated with this foliation are the A'^ — 3 functionally independent cross ratios of the 
oscillator phases. No further reduction is possible, in general; numerical experiments on models 
of Josephson junction arrays suggest that the invariant manifolds often contain three-dimensional 
regions of neutrally stable chaos. 



Large arrays of coupled limit-cycle oscillators 
have been used to model diverse systems in 
physics, biology, chemistry, engineering and so- 
cial science. The special case of phase oscillators 
coupled all-to-all through sinusoidal interactions 
has attracted mathematical interest because of its 
analytical tractability. About 20 years ago, nu- 
merical experiments revealed that these systems 
display an exceptionally simple form of collective 
behavior: for all N >3, where TV is the number of 
oscillators, cdl trajectories are confined to meuii- 
folds with N — 3 fewer dimensions than the state 
space itself. Several insights have been obtained 
over the past two decades, but it has remained an 
open problem to pinpoint the symmetry or other 
structure that causes this non-generic behavior. 
Here we show that group theory provides the ex- 
planation: the governing equations for these sys- 
tems arise naturally from the action of the group 
of conformal mappings of the unit disc to itself. 
This link unifies and explains the previous numer- 
ical cind cinalytical results, cind yields new con- 
stants of motion for this class of dynamical sys- 
tems. 



I. INTRODUCTION 

When a nonlinear system shows unexpectedly simple 
behavior, it may be a clue that some hidden structure 
awaits discovery. 

For example, recall the classic detective story [1] that 
began in the 1950s with the work of Fermi, Pasta, and 
Ulam [2, 3, 4]. In their numerical simulations of a chain 
of anharmonic oscillators, Fermi et al. were surprised 
to find the chain returning almost perfectly, again and 



'Electronic address: sam255@cornell.edu 



again, to its initial state. The struggle to understand 
these recurrences led Zabusky and Kruskal [5] to the 
discovery of solitons in the Korteweg-deVries equation, 
which in turn sparked a series of results showing that 
this equation possessed many conserved quantities in 
fact, infinitely many [6]. Then several other equations 
turned out to have the same properties. At the time 
these results seemed almost miraculous. But by the mid- 
1970s the hidden structure responsible for all of them — 
the complete integrability of certain infinite-dimensional 
Hamiltonian systems [7] — had been made manifest by the 
inverse scattering transform [8, 9] and Lax pairs [10]. 

Something similar, though far less profound, has been 
happening again in nonlinear science. The broad topic is 
still coupled oscillators, but unlike the conservative oscil- 
lators studied by Fermi et al., the oscillators in ques- 
tion now are dissipative and have stable limit cycles. 
This latest story began around 1990, when a few re- 
searchers noticed an enormous amount of neutral stabil- 
ity and seemingly low-dimensional behavior in their sim- 
ulations of Josephson junction arrays — specifically, ar- 
rays of identical; ovcrdampcd junctions arranged in series 
and coupled through a common load [11, 12, 13, 14, 15]. 
Then, just a year ago, Antonsen et al. [16] uncovered 
similarly low-dimensional dynamics in the periodically 
forced version of the Kuramoto model of biological os- 
cillators [17, 18, 19]. This was particularly surprising 
because the oscillators in that model are non-identical. 

As in the soliton story, these numerical observations 
then inspired a series of theoretical advances. These in- 
cluded the discovery of constants of motion [20, 21], and 
of a pair of transformations that established the low- 
dimensionality of the dynamics [20, 21, 22, 23, 24, 25]. 
But what remained to be found was the final piece, the 
identification of the hidden structure. Without it, it was 
unclear why the transformations and constants of motion 
should exist in the first place. 

In this paper we show that the group of Mobius trans- 
formations is the key to understanding this class of dy- 
namical systems. Our analysis unifies the previous treat- 
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mcnts of Josephson arrays and the Kuramoto model, and 
clarifies the geometric and algebraic structures responsi- 
ble for their low-dimensional behavior. One spin-off of 
our approach is a new set of constants of motion; these 
generalize the constants found previously, and hold for a 
wider class of oscillator arrays. 

The paper is organized as follows. To keep the treat- 
ment self-contained and to estabUsh notation, Section II 
reviews the relevant background about coupled oscilla- 
tors and the Mobius group. In Section III we show how 
to use Mobius transformations to reduce the dynamics of 
oscillator arrays with global sinusoidal coupling, a class 
that includes the Josephson and Kuramoto models as 
special cases. The reduced flow lives on a set of invariant 
three-dimensional manifolds, arising naturally as the so- 
called group orbits of the Mobius group. The results ob- 
tained in this way are then compared to previous findings 
(Section IV) and used to generate new constants of mo- 
tion via the classical cross ratio construction (Section V). 
We explore the dynamics on the invariant manifolds in 
Section VI, and show that the phase portraits for resis- 
tively coupled Josephson arrays are filled with chaos and 
island chains, reminiscent of the pictures encountered in 
Hamiltonian chaos and KAM theory. 



II. BACKGROUND 
A. Reducible systems with sinusoidal coupling 



for j = I, . . . ,N. By direct substitution, one can check 
that the resulting functions (t>j{t) simultaneously sat- 
isfy all equations in (1) as long as the three vari- 
ables $(t),7(t) and Q{t) satisfy a certain closed set of 
ODEs [21]. 

Watanabe and Strogatz also noted that the same trans- 
formation can be used to reduce any system of the form 

4>,^ fe'^^ +g + fe-''^^ (3) 

for j = 1, . . . ,N, where / is any smooth, complex- valued, 
27r-periodic function of the phases (pi, . . . (Here the 
overbar denotes complex conjugate. Also, note that g 
has to be real- valued since is real.) The functions / 
and g are allowed to depend on time and on any other 
auxihary state variables in the system, for example, the 
charge on a load capacitor or the current through a load 
resistor for certain Josephson junc;tion arrays. The key 
is that / and g must be the same for all oscillators, and 
thus do not depend on the index j. We call such systems 
sinusoidally coupled because the dependence on j occurs 
solely through the first harmonics e*"^^ and e"*"^^ . 

Soon after the transformation (2) was reported, 
Goebel [22] observed that it could be related to fractional 
linear transformations, and he used this fact to simplify 
some of the calculations in Ref. [21]. At that point, re- 
search on the reducibility of Josephson arrays paused for 
more than a decade. The question of why this particu- 
lar class of dynamical systems (3) should be reducible by 
fractional linear transformations was not pursued at that 
time, but will be addressed in Section III. 



The theory developed here was originally motivated by 
simulations of the governing equations for a series array 
of N identical, overdamped Josephson junctions driven 
by a constant current and coupled through a resistive 
load. As shown in Tsang et al. [11], the dimensionless 
circuit equations for this system can be written as 



n-{b+l) cos ( 



1 ^ 



COS (pk 



(1) 



for j = I, . . . , N . The physical interpretation need not 
concern us here; the important point for our purposes is 
that this set of A'^ ordinary differential equations (ODEs) 
displayed low-dimensional dynamics. The same sort of 
low-dimensional behavior was later found in other kinds 
of oscillator arrays [14] as well as in Josephson arrays 
with other kinds of loads [12, 13, 15]. 

Building on contributions from several teams of re- 
searchers [11, 12, 13, 14, 15], Watanabe and Strogatz [21] 
showed that the system (1) could be reduced from A'' 
ODEs to three ODEs, in the following sense. Consider 
a time-dependent transformation from a set of constant 
angles 6j to a set of functions ^j{t), defined via 



tan 



1 - 7(i) 



tan 



9j - e{t) 



(2) 



B. Ott-Antonsen ansatz 

Ott and Antonsen [23, 25] recently reopened the issue 
of low-dimensional dynamics, with their discovery of an 
ansatz that collapses the infinite-dimensional Kuramoto 
model to a two-cliincnsional system of ODEs. 

To illustrate their ansatz in its simplest form, let us 
apply it to the class of identical oscillators governed by 
Eq.(3), in the limit N ^ oo. (Note that this step in- 
volves two simplifying assumptions, namely, that N is 
infinitely large and that the oscillators are identical. The 
Ott-Antonsen ansatz applies more generally to systems of 
non-identical oscillators with frequencies chosen at ran- 
dom from a prescribed probability distribution — indeed, 
this generalization was one of Ott and Antonsen's major 
advances — but it is not needed for the issues that we wish 
to address.) In the limit N ^ oo, the evolution of the 
system (3) is given by the continuity equation 

dp djpv) 

dt d(f> 



= 



(4) 



where the phase density p{4',t) is defined such that 
p{(l),t)d(p gives the fraction of phases that lie between 
(j) and cf) + dd) at time t, and where the velocity field is 
the Eulerian version of (3): 



v{cP,t) = fe'^ + g + fe-"^. 



(5) 
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Our earlier assumptions about the coefficient functions / 
and g now take the form that / and g may depend on t 
but not on 0. The time-dependence of / and g can arise 
either explicitly (through external forcing, say) or implic- 
itly (through the time-dependence of the harmonics of p 
or any auxiliary state variables in the system). 

Following Ott and Antonsen [23], suppose p is of the 
form 



27r\ 



^(a(t)"e*"* + a(t)"e-*"'^) 



(6) 



for some unknown function a that is independent of </>. 
(Our definition of a is, however, slightly different from 
that in Ott and Antonsen [23]; our a is their a.) Note 
that (6) is just an algebraic rearrangement of the usual 
form for the Poisson kernel: 



l-r2 



27r 1 - 2r cos(^ - $) + 
where r and $ are defined via 



a = re 



(7) 



(8) 



In geometrical terms, the ansatz (6) defines a submani- 
fold in the infinite-dimensional space of density functions 
p. This Poisson submanifold is two-dimensional and is 
parametrized by the complex number a, or equivalently, 
by the polar coordinates r and 

The intriguing fact discovered by Ott and Antonsen is 
that the Poisson submanifold is invariant: if the density 
is initially a Poisson kernel, it remains a Poisson kernel 
for all time. To verify this, we substitute the velocity field 
(5) and the ansatz (6) into the continuity equation (4), 
and find that the amplitude equations for each harmonic 
gjn<^ are simultaneously satisfied if and only if a{t) evolves 
according to 



a = i{f + ga + fa^). 



(9) 



This equation can be recast in a more physically mean- 
ingful form in terms of the complex order parameter, de- 
noted by {z) and defined as the centroid of the phases <j) 
regarded as points e"^ on the unit circle: 



^2fr 

Jo ' 



e**p(<^,i)d<^. 



(10) 



By substituting (6) into (10) we find that (z) = a for all 
states on the Poisson submanifold. Hence, {z) satisfies 
the Riccati equation 



{z)=i{f + g{z) + f{zr). 



(11) 



When / and g are functions of (z) alone, as in 
mean- field models, Eq.(ll) constitutes a closed two- 
dimensional system for the flow on the Poisson submani- 
fold. More generally, the system will be closed whenever 
/ and g depend on p only through its Fourier coefficients. 



We will show this explicitly in Subsection VB, by finding 
formulas for all the higher Fourier coefficients in terms of 
a, and hence in terms of (z). (However, as we will see, 
things become more complicated for states lying off the 
Poisson submanifold. Then (z) no longer coincides with 
a and the closed system becomes three dimensional, in- 
volving tp as well as a.) 

The work of Ott and Antonsen [23] raises several ques- 
tions. Why should the set of Poisson kernels be invari- 
ant? What is the relationship, if any, between the ansatz 
(6) and the transformation (2) studied earlier? Why 
does (2) reduce equations of the form (3) to a three- 
dimensional flow, whereas (6) reduces them to a two- 
dimensional flow? 

As we shall see, the answers have to do with the prop- 
erties of the group of conformal mappings of the unit disc 
to itself. Before showing how this group arises naturally 
in the dynamics of sinusoidally coupled oscillators, let us 
recall some of its relevant properties. 



C. Mobius group 

Consider the set of all fractional linear transformations 
: C ^ C of the form 



F{z) 



az - 



CZ + d 



(12) 



where a,h,c and d are complex numbers, and the nu- 
merator is not a multiple of the denominator (that is, 
ad — be ^ Q). This family of functions carries the struc- 
ture of a group. The group operation is composition of 
functions, the identity element is the identity map, and 
inverses are given by inverse functions. 

Of most importance to us is a subgroup G — which we 
refer to as the Mobius group — consisting of those frac- 
tional linear transformations that map the open unit 
disc n = {z & C : \z\ < 1} onto itself in a one-to- 
one way. These transformations and their inverses are 
analytic on D and map its boundary (the unit circle 
S''^ = {z € C : \z\ = 1}) to itself. All such automor- 
phisms of the disc can be written [26] in the form 



F{z) = 



1 — ' 

1 — az 



(13) 



for some (/? € M and a G D. The Mobius group G is in 
fact a three-dimensional Lie group, with real parameters 
ip, Re(a), and Im(Q!). 

However, it turns out that a different parametrization 
of G will be more notationally convenient in what follows, 
in the sense that it simplifies comparisons between our 
results and those in the prior literature. Specifically, we 
will view a typical element of G as a mapping M from 
the unit disc in the complex w-plane to the unit disc in 
the complex ^-plane, with parametrization given by 



M{w) = - 

1 



(14) 
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where a e D and V € M. Note that the inverse mapping Prom (15), we get 



w 



■ az 



(15) 



has an appearance closer to that of the standard 

paramctrization (13). 

A word about terminology: our definition of the 
Mobius group is not the conventional one. Usually this 
term d(mot(^s the larger group of all fractional linear 
transformations (or bilinear transformations, or linear 
fractional transformations), whereas we reserve the ad- 
jective Mobius for the subgroup G and its elements. 
Thus, from now on, when we say Mobius transforma- 
tion we specifically mean an clement of the subgroup G 
consisting of analytic automorphisms of the unit disc. 



III. MOBIUS GROUP REDUCTION 

In this section we show that if the equations for the 
oscillator array are of the form (3), then the oscillators' 
phases {t) evolve according to the action of the Mobius 
group on the complex unit circle: 



(16) 



for j = 1,...,A'", where is a one-parameter family 
of Mobius transformations and 6j is a constant (time- 
independent) angle. In other words, the time-t flow map 
for the system is always a Mobius map. 

Incidentally, this result is consistent with a basic topo- 
logical fact: we know that different oscillators can- 
not pass through each other on under the flow, so 
we expect the time-t flow map to be an orientation- 
preserving homeomorphism of onto itself — and indeed 
any Mobius map is. 

We begin the analysis with an algebraic method similar 
to that in Goebel [22]. Then, in Sections IIIB and IIIC, 
we adopt a geometrical perspective and show that it an- 
swers several questions left open by the first method. 



A. Algebraic Method 

Parametrize the one-parameter family of Mobius trans- 
formations as 



Mt{w) 



1 -|- ae'^'^w 
where \a{t)\ < 1 and ip{t) e K, and let 



(17) 



(18) 

To verify that (17) gives an exact solution of (3) — subject 

to the constraint that the Mobius parameters a(t) and 
'ijj{t) obey appropriate ODEs, to be determined — we com- 
pute the time-derivative of (f>j{t) = — ilogMt(wj), keep- 
ing in mind that Wj is constant: 

• ipe^'^Wj — ia {ia — aip)e^'^Wj 



■' 1 - ae'»'j 
which when substituted into (19) yields 

ip + ioLOL — aiioL — aip) 



1 — \oi\ 



(20) 



Re-"^^ (21) 



where R ~ {ia — at/j) /(I — |ap). 

Note that Eq.(21) falls precisely into the algebraic form 
required by (3). Thus, to derive the desired ODEs for 
a{t) and ijj{t), we now subtract (21) from (3) to obtain 
N equations of the form = Cie"^^ + Cq + C.ie-*'^^ 
fov j = 1, . . . ,N. If the system contains at least three 
distinct oscillator phases, then Ci, Co, and C_i must 
generically be zero. Explicitly, 



/ = 



a- 



tjj + iaa — a{ia — aip) , , 
5= ^ VTvi • (22) 



The system (22) can be algebraically rearranged to give 

a = i{fa^ + ga + f) (23a) 
i^ = fa + g + fa. (23b) 

Equations (23a) and (23b) have been derived previ- 
ously; they appear as Eqs.(lO) and (11), respectively, in 
Pikovsky and Rosenblum's work [24], where they were 
derived by applying the transformation (2). Both their 
approach and the one above are certainly quick and clean, 
but they require us to guess the transformation ahead 
of time, and reveal little about why this transformation 
works. 

Incidentally, observe that under the change of variables 

Zj = e^'^' , (3) becomes 



= i{fZj+gzj+f). 



(24) 



1 -|- ae^^Wj 



Equation (24) is a Riccati equation with the form of 
(23a) — another coincidence that seems a bit surprising 
when approached this way. In the following subsection, 
we will see how these Ricatti equations emerge naturally 
from the infinitesimal generators of the Mobius group. 



B. Geometric Method of Finding a 

Now we change our view of Mobius maps slightly. In- 
stead of thinking of M as a map from the ui-plane to the 
^;-plane, we view it as a map from the 2;-plane to itself. 
This requires a small and temporary change in notation, 
but it makes things clearer, especially when we start to 
discuss differential equations on the complex plane. 

We begin by recalling some basic facts and defini- 
tions. Suppose the coupled oscillator system contains 
just three distinct phases among its oscillators. Then 
by a property of Mobius transformations, there exists 
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a unique Mobius transformation from any point Z\ 



(e 



i6\ „i02 pi8 



e^^) to any other point Z2 = (e 



101 



2*02 ^g*03) in 

the state space 5^ x x S^. If the system instead con- 
tains only one or two distinct phases, many Mobius trans- 
formations will take Zi to Z2, so we can still reach every 
point of the phase space from every other point. However, 
if the system contains more than three distinct phases, 
say N, then there is not in general a Mobius transfor- 
mation that transforms Zi = (e*^% e*^^, e*^'', . . . , e'^") to 
Z2 = (e''^i , e"^2^ e*"^S . . . , e''^"); only some points are ac- 
cessible from Zi, while others are not. 

In the language of group theory, we say that Z2 is in 
the group orbit of Zi if there exists a Mobius map M such 
that Z2 = M(zi). Then, as a direct consequence of the 
fact that Mobius maps form a three-parameter group G 
under composition, the group orbits of G partition the 
phase space into three-dimensional manifolds (when the 
phase space is at least three-dimensional). 

To compute infinitesimal generators for G, we compute 
the time derivatives of the three one-parameter families 
of curves corresponding to the three parameters of G: 
ip, Re(a) and Im(a). Each of the three families is ob- 
tained from the Mobius transformation by setting two of 
the three parameters to zero, and leaving the remain- 
ing parameter free. For example, if we set t = at 
z = {zi, . . . , zn), these three families are 



Mi(z) = e'*z 
M2(z)= 

Msiz) = 



1 - tz 

z + it 



(25) 



itz 



where Mi{z) is written in place of {Mi{zi), . . . , A'/i(zjv)) 
and likewise for M2(z) and Ms{z). We continue using 
this shorthand in subsequent equations, writing h{z) in 
place of {h{zi), . . . ,h{zN)) for any one-parameter func- 
tion h. 

The time derivatives of the curves in (25) evaluated at 
t = then give a set of infinitesimal generators for G: 



Vi = iz 

V2 = z'^-1 

V3 = iz^ + i. 



(26) 



V = ciVi + • • • + CnVn whcrc the coefficients depend 
only on time t, then the trajectory of the dynamics z = v 
from any initial condition Zq can be expressed in the form 
{At{zo)} for a unique family {At} C L parameterized 
by t. Since the Mobius group is a complex Lie group, 
this result can be applied directly to conclude (27) has 
the solution z{t) = Mt(zo) where {M^} is a unique one- 
parameter family of Mobius transformations. 

Although we have so far assumed that the components 
Zk of z lie on the complex unit circle, both (17) and 
(27) extend naturally to all of C^. This implies that 
Z() = must evolve as z{t) = Mi{0) for some family 
{Ml}. However, Eq. (17) shows that M(0) = a for all 
M e G. So z{t) = Mt{0) = a for all t, meaning that 
a{t) satisfies (27). Since (27) is just a rewriting of (24), 
the dynamics (23a) for a that we derived earlier are now 
placed in a geometrical context. This approach reveals 
that a{t) is just the image of the origin under a one- 
parameter family of Mobius maps, applied to any one 
complex plane of C"'^. 

It is even more illuminating to compute the infinitesi- 
mal generators within the A''-fold torus of phase val- 
ues, i.e., the quantities Uk = — logMfe(e"^)|t=o. These 
turn out to be 



Ml = (!,..., 1) 
U2 = 2 sin 
U3 = 2 cos <p. 



(28) 



When expressed in terms of these infinitesimal genera- 
tors, the equation of motion (27) becomes 



4> = g + {2 sin 4>)hi + (2 cos 0)/i2 



(29) 



which is precisely what we earlier referred to as a si- 
nusoidally coupled system (3), and whose solution must 
therefore be of the form <pt = — ilogMt(e'^) for some 

Mt e G. 

This calculation finally clarifies what is so special 
about sinusoidally coupled systems (3) : they are induced 
naturally by a fiow on the Mobius group. This fact un- 
derlies their reducibility and all their other beautiful (but 
non-generic) properties. 



C. Geometric Method of Finding tp 



Note that these three generators point out into the full 
AT-dimensional complex space C^, as expected. 

Meanwhile, if we substitute / = —ihi + /12 (where 
hi and /12 are real functions) into the original Riccati 
dynamics (24), we can rewrite this equation of motion in 
terms of the three infinitesimal generators: 



z = izg + (z^ — l)hi + {iz"^ + i)h2- 



(27) 



The implication of the rewritten form (27) is then given 
by a theorem from Lie theory: if L is a Lie group acting 

on a submanifold with linearly independent infinitesimal 
generators V\,. . . ,Vn, and v is a vector field of the form 



We turn next to the dynamics of V'- As we will show 

in the next section, the action of the Mobius transforma- 
tion involves a clockwise rotation of the oscillator phase 
density p{(j), t) by arg(Q!) — i/i and a counterclockwise rota- 
tion by arg(a). Hence, ipit) may be viewed as the overall 
counterclockwise rotation of the distribution at time t 
relative to the initial distribution at t = 0. 

To support this interpretation, we show here that ip 
equals the average value of the vector field on the circle, 
given by 



(30) 
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Observe the right side of the integrand (19) has two 
terms; 



Riiw) = 

R2{w) = 
By Cauchy's formula, 



ipe^'^w — ia 
{ia 



1 + ae»^ 



w 



27ri J SI 
So {(j)) simpUfies to 

(0) = 



1 



dw 



w 



i?2(0) = 0. 



Ri{w) 



dw 



w 



(31) 



(32) 



(33) 



Note that i?i(w) has a pole in the unit disc, so we make 
the change of variables w w^^ to move this pole out- 
side the circle. Evaluating the resulting integral yields 



1 /" „ . s dw 
— / RAw) — 
zmJsi w 



1 f — icee 
2^1 J gi 1 + ae 



— iae ^^w dw 



-ill) 



w w 



^ (34) 



which completes the demonstration that (0) = ijj. 

We can now go back and evaluate the average vec- 
tor field in a different way to find the differential equa- 
tion that governs ip{t). Differentiating <p = —ilogMt{w) 
with respect to time and substituting the result into 
ip = ^ Jgi ^d6, we obtain 



liri J SI 



Mt{w) dw 
2'Ki J SI Mt{w) iw ' 



(35) 



Since obeys the Ricatti equation, we can eliminate 

Mt in the numerator above to get 

1 /" d 

i' = wz: I ifMt{w) + g + fMt{w)-^) — . (36) 



w 



There are three integrals to evaluate here. The third one 
involves a term Mt{w)~^ which has a pole inside the unit 
circle, so we do the same change of variables as before, 
w w~^, to move the pole outside. The corresponding 

integral then simplifies to 

_i dw 



2TTi J SI 



Mt{w)- 



w 



-I - 

2TTi J SI 1 



dw 



+ ae 



w w 



a (37) 



where the final integration follows from Cauchy's for- 
mula. Similiarly, we use Cauchy's formula to integrate 
the first and second terms of the integrand in (36), and 
thereby obtain the desired differential equation for ijj, 
thus rederiving (23b) found earlier. 



IV. CONNECTIONS TO PREVIOUS RESULTS 



Relation to the Watanabe-Strogatz 
Transformation 



It is natural to ask how the trigonometric transforma- 
tion (2) used in earlier studies [20, 21, 24] relates to the 



Mobius transformation (17) used above. As we will see, 
(2) may be viewed as a restriction of (17) to the complex 
unit circle. 

First, by trigonometric identities, we have 



tan 



1 



'l + e»W-*)' 



(38) 



To connect this to Mobius transformations, consider 
what happens when we apply the map defined by (17) 
to a point w = e'^ on the unit circle. Since the image 
is also a point on the unit circle, it can be written as 



M(e*^) = e**^ for some angle </>. Next let a 
divide both sides of (17) by e'*. Thus 



re 



and 



i(6>-0) 



+ r 



1 -I- j.gi{e-e) 



(39) 



where Q = ^ — tjj. Substitution of (39) into the right side 
of (38) gives 



tan 





1 — r / 


2 


~ l + r\ 



1 



=*(e-e) 



l + e»(^-e) 



(40) 



By the identity (38), Eq.(40) is equivalent to (2) with 
7 = _2r/(l-|-r2). 

We can now see how the Mobius parameters a and 
operate on the set of e*^ in C. From the relationships 
between 8, 7, $ and the Mobius parameters, the ini- 
tial phase density is first rotated clockwise around by 
arg(a) — ip, then squeezed toward one side of the cir- 
cle as a function of \a\, and afterwards rotated counter- 
clockwise by arg(a). The squeeze, which takes uniform 
distributions to Poisson kernels, can be thought of as a 
composition of inversions, dilations and translations in 
the complex plane. 



B. Invariant Manifold of Poisson Kernels 

In Section II B and in a previous paper [27] , we used the 
Ott-Antonsen ansatz (6) to show that systems of identi- 
cal oscillators with global sinusoidal coupling contain a 
degenerate two-dimensional manifold among the three- 
dimensional leaves of their phase space foliation. This 
two-dimensional manifold, which we called the Poisson 
submanifold, consists of phase densities p{(f), t) that have 
the form of a Poisson kernel. We now rederive these re- 
sults within the framework of Mobius transformations. 

Let T denote one instance of the transformation (2); 
in other words, fix the parameters 7 and 6 and let 
4> = T{9). Let fj, denote the normalized uniform measure 
on S^; thus 



di^{e) 



1 

2^" 



-dO. 



(41) 



The transformation T maps /i to the measure T^/i, and, 
by the usual formula for transformation of single- variable 
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measures, we have d{T^,fi){(f>) = i^T^^ {(/>)' d(j), where the 
prime denotes differentiation by (f>. From this equation it 
follows that d{T^ij){(t>) has the form of the Poisson kernel, 
because the inverse of the Mobius transformation (17) is 



M-'^{z)=e 



— p-ii' _ 



(42) 



1 — az 
which implies 

r-i((/)) = -il)- ilog(e"^ - a) + i log(l - ae*"^). (43) 

Then by differentiation and algebraic rearrangement, we 
obtain 



T-\<i>)' 



1 - 2rcos((/) - <I>) + r2 



(44) 



The integral of r-i(0)' over [0,27r) is 27r, so d{T^ii){(l)) 
is indeed a normalized Poisson kernel. 

Finally, if the phase distribution d{T^, fi) {(p) / d(j) ever 
takes the form of a Poisson kernel with parameters r = tq 
and $ = i>o, then we can set r(0) = ro, <I>(0) = i>o 
and dy,{9) — ^dO, and the above calculation shows that 
d{T^lj){4') / dfj) remains a Poisson kernel for all future and 
past times. Hence, the set of normalized Poisson ker- 
nels constitutes an invariant submanifold of the infinite- 
dimensional phase space. 

The above demonstration also reveals that the Poisson 
submanifold has dimension k + 2 where k is the number 
of state variables besides a, and the oscillator phases. 
More concretely, it implies that when the system lies on 
the Poisson submanifold, we can write a as depending 
only on a; it is not possible to require a to depend on ip 
in any real coupling scheme. 

To see this, we first consider the case in which the sys- 
tem is closed and there are no additional state variables. 
Suppose a does not depend only on a. Then some of the 
state space trajectories cross when projected onto the 
unit disc of a values. At the point of any crossing, the 
phase density p{(t>,t) has multiple d values. But by (44), 
the phase density depends only on a, so there is nothing 
in the state space that can distinguish between the differ- 
ent a values at that point. Hence, a must be expressible 
in terms of a alone. By an analogous argument, a is also 
independent of on the Poisson submanifold when there 
are k other state variables besides the oscillator phases 
and Mobius parameters. 

On the other hand, if the time-dependence of a arises 
only via a dependence on a, then r and $ decouple 
from ij) and the dynamics are two-dimensional regard- 
less of whether the system is evolving on the Poisson 
submanifold or not. Observe that we can always force V'- 
indepencleiicie for a by throwing away enough information 
about the locations of the other phases. For instance, in 
the extreme, we may simply make / and g constant. 

Finally, even when a does not depend solely on a, the 
dynamics still may be two dimensional. For example, in 
the case of completely integrable systems [20], the vari- 
ables r and <I> — ■0 decouple from $ to foliate the phase 
space with two-dimensional tori. 



V. CHARACTERISTICS OF THE MOTION 

A. Cross Ratios as Constants of Motion 

The reduction of (3) by the three-parameter Mobius 
group suggests that the corresponding system of coupled 
oscillators should have — 3 constants of motion. As we 
will see, these conserved quantities are given by the cross 
ratios of the points Zj = e"^^ on S^. Recall from complex 
analysis [28] that the cross ratio of four distinct points 
zi,Z2, z^, 2:4 e C U {00} is 



{zi,Z2, Z3, Zi) 



Zl - Z3 Z2 - Z4 
Z\ — Z4 Z2 — Z3 



(45) 



This quantity is conserved under Mobius transforma- 
tions: for all a and V, {M{zi),M{z2),M{z3),M{z4)) = 
{zi, Z2, Z3, Z4). Hence, the A^!/(A — 4)! cross ratios of 
the A^ oscillator phases remain constant along the tra- 
jectories in phase space. We denote the constant value 
of (-21, 2^2, -23, ^4) as A1234. Of course, we could have de- 
fined the cross ratio for four-tuples of non-distinct points 
as well, but these quantities are trivially conserved re- 
gardless of the dynamics and hence do not reduce the 
dimension of the phase space. 

To show that exactly A^ — 3 of the cross ratios 
are independent, consider the sequence: {{zi, Z2, Z3, Z4), 
{Z2, Z3, -24, -25), . . . , (-2JV-3, -2Ar-2, zn-1, -2jv)}- Each cross 
ratio in the sequence includes a new point not in the 
cross ratios preceding it and therefore must be indepen- 
dent of them. Hence, there are at least A^— 3 independent 
cross ratios. With a bit more work (see the Appendix), 
we can also confirm that the rest of the cross ratios are 
functionally dependent on these A^ — 3 integrals. 

Since the state space of the phases is an A'-fold torus 
of real variables, we expect that each of the constants 
of motion can be expressed in terms of real functions 
and variables. Indeed, if zi, Z2, -23, -24 lie on the unit 
circle, then the cross ratio (zi, Z2, Z3, Z4) lies on MU {00}. 
We see this explicitly by pulling out e5('*i+'^-'^^ from the 
factor {e''l'' - e"''^ ) of (e'*i , e*<^^ e"^^ e*'P^), and likewise 
for the other three factors, and then canceling the factors 
g^i{4>i+<t'2+4>3+<t'4) [j^ ^]^Q numerator and denominator to 
find 



(g«<)!.l^gi<^2^gi<A3^gi<?i4) 



S13S: 



13'J24 



S-iaS- 



14'-'23 



where 



Sij = sin 



(46) 



(47) 



This way of writing the cross ratio also suggests a rela- 
tionship with the constants of motion reported by Watan- 
abe and Strogatz [20, 21] for completely integrable sys- 
tems (those with / = ^e'^^{z) and g = 0, where (-2) is the 
phase centroid (10)). These constants of motion, which 
we will call WS integrals, take the form 



<S'l2'S'23 ■ ■ ■ 



(48) 
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where any permutation of the indices generates another 
WS integral. As previously demonstrated [21], exactly 
A'' — 2 of the N\ index permutations of (48) are function- 
ally independent. 

As we might anticipate, the WS integrals imply that 
the cross ratios are constants of motion: consider two dis- 
tinct WS integrals / = SikSkiSijH and /' = SuSikSkjH, 
where 11 denotes the remaining product of factors. As- 
sume n is the same for both / and /'. Then ///' = 
—\jki- Since «, j, fc, / are arbitrary, we can generate all 
cross ratios via this procedure. 

Additionally, if a single WS integral holds for a system 
in which the cross ratios are invariant, then all WS inte- 
grals hold, since we can arbitrarily permute the indices 
of the first WS integral by sequences of transpositions of 
the form I = —Xijkil' in which I and k are interchanged. 



Now, M(w)™(l + P{w)) is analytic on the open disc D 
and M(0)™(1 -|- -P(O)) = a™. Meanwhile, the remaining 
term of the integrand of (52) has the complex conjugate 



M{w) P{w) 
w 



1 + ae'^^w 
e^^w + a 



'P(w) 
w 



(53) 



which features an order- 1 pole act w = and an order-m 

pole at w = —e^'^^'a. The first residue evaluates to zero, 
while the second is given by 



(to - 1)! dw™- 



(1 + ae'l'w) 



w 



(54) 



Therefore, {z"^) is equal to a™ added to the complex 
conjugate of this second residue: 



B. Fourier Coefficients of the Phase Distribution i^"^) — o;™ + ^ ^ 



2\k+l 



When we introduced / and g in Section II, we required 
that they depend on the phases only through the Fourier 

coefficients of the phase density p{(f>,t). Since the cen- 
troid (10) is the Fourier coefficient corresponding to the 
first harmonic e~"^, this condition is met by standard 
Kuramoto models, Josephson junction series arrays, laser 
arrays and many other well-studied systems of globally 
coupled oscillators. 

Our goal now is to show that this condition implies the 
closure of (23), in the sense that a and ip depend only 
on a and ijj. To do so, we will show that the Fourier 
coefficient of all higher harmonics e"""'^ for any integer 
TO may be expressed in terms of a and t/j- 

For a fixed measure fi(6) on [0, 27r) and a transforma- 
tion T{9) = —i log M{0) of this measure via the Mobius 
map M, the Fourier coefficient of e"""''* is given by 



{zn = [ 



(r,M)(0) = / M{e^y^di^{0). (49) 



We use the notation {z™) as a reminder that {z) is the 
phase centroid. 

We assume that we can take a Fourier expansion of 
IJ.{e), so 



dn{e) 



1 

2^ 



E 



(50) 



where the constants c„ are independent of 9. Since the 
phase distribution must be real and normalized, we know 
that c_„ = c„ and cq = 1, so we can write 



d,= -[l + P{v:) 



/ w 



(51) 



where w — e*^ and P{w) = X^^i CnW^. The formula for 
{z"^) then becomes: 

lz^) = —[ M{wr(l + P(w) + F^^—. (52) 
2Tn Js^ \ J w 



fe=o 



E(-i) 



k\ 



n=0 



c„+fe+ie'('"+")'^a". (55) 



For example, the centroid may be written in terms of a 
and ij} as 

oo 

{z)=a + (|a|2 - 1) 5^(-l)"c„e'"'^a"-^ (56) 

This calculation reveals what is so special about the 
Poisson submanifold. Recall from Section IV B that Pois- 

son kernels arise when we take /i to be the uniform mea- 
sure. Then c„ = for all n ^ and (z) — a. In this 
exceptional case, the centroid simply evolves according 
to the Riccati equation (11) and the dynamics of a and 
■0 decouple in Eqs. (23a), (23b). (A similar observation 
about the crucial role of the uniform measure here was 
made by Pikovsky and Rosenblum [24]. The centroid 
evolution equation (23a) on the Poisson submanifold was 
first written down by Ott and Antonsen; see Eq.(6) in 
Ref. [23].) 

But for the generic case of states lying off the Poisson 

submanifold, (z) is no longer equal to a and the reduced 
dynamics become fully three-dimensional, due to the cou- 
pling between a and tl) induced by the relation (56) and 
the dependence of / and g on (z) and the higher Fourier 
coefficients. In the next section we will explore some of 
the possibilities for such three-dimensional flows. 



VI. CHAOS IN JOSEPHSON ARRAYS 

Although the leaves of the foliation imposed by the 

Mobius group action are only three-dimensional, they of- 
ten contain chaos for commonly studied / and g [14, 21]. 
In this section, we showcase this phenomenon by spe- 
cializing to the case of a resistively-loaded series array of 
overdamped Josephson junctions. 
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FIG. 1: The qualitative trend of chaos observed in the first 
quadrant of the b-fl parameter plane is indicated by the 
shaded gradient. As the shade darkens near the bifurcation 
curve Q = b, chaos fills increasingly larger regions of the 
submanifolds containing the sinusoidal initial distributions. 
Points (A) and (B) are chosen as (1/20, 3/4) and (17/10, 1), 
respectively. Representative Poincare sections for these points 
are shown in Fig. 2 and Fig. 3. The region 6 < is grayed 
out to represent that negative values of b are not physical. 

In several previous studies of sinusoidally coupled os- 
cillators in the continuum limit, it was found that under 
certain conditions, the Fourier harmonics of the phase 
density p{(j), t) evolved as if they were decoupled, at least 
near certain points in state space [14, 29, 30]. In the 
spirit of these observations, we can get a sense for how 
individual harmonics contribute to the chaos by start- 
ing the system (23) from sinusoidal phase densities with 
different wavenumbers n. 

To be more precise, we choose an initial density 

p(0,O) = ^(l+cosn0). (57) 

At t = 0, we choose a = V' = so that Mt in Eq. (17) is 
simply the identity map, and the time-dependent change 
of variables e"^ = Mt(e*^) reduces to (f) — 9^ initially. 
Thus, the corresponding density of 9 is 

a^{9) ^ ^ {I + cos n9). (58) 
27r 

This density is independent of time, just as the angles 9j 
were in the finite- TV case. 

Next we flow the density forward by e"^ — Mt(e*^), 
where the Mobius parameters a(t)^'ip(t) satisfy the re- 
duced flow (23). Then, by our earlier results, the re- 



sulting density p{(j), t) automatically satisfies the govern- 
ing equations (4), (5). The three-dimensional plot of 
Re{a{t)), Im(a(t)) and ^p{t) indicates how such a single- 
harmonic density evolves in time, revealing for example 
whether it exhibits chaos, follows a periodic orbit, or ap- 
proaches a fixed point. 

To ease the notation, from now on we write a in Carte- 
sian coordinates as 

a — X + iy. (59) 

Then the reduced flow (23) becomes 

X — —uy + Im(/)(1 — — j/^) 

y = ux + Re{f){\-x'^ -y"^) (60) 

ij} — u 

where 

u^2xRe{f) + g~2ylui{f). (61) 

We immediately see that for every fixed point of this 
system, jaj = 1 and ip is arbitrary. If for some change 
of state variables (^{x,y,ip)y 'yls;,^, and ^{x,y,'ip), the 
ODEs C and 77 constitute a closed two-dimensional system 
and ^ receives all of its i-dependence through C and 77, 
then there could be other fixed points for the physical 
system, namely where C = ^ = but =/= 0- Examples 
of the second type of fixed point include the splay states 
found on the Poisson submanifold [11, 30]. 

As discussed in Section II A, series arrays of Josephson 
junctions with a resistive load have dynamics given by 
Eqs. (1), (4), and (5), with / -{b+l)/2 and g = fl + 
Re(z), where b and f2 are dimensionless combinations of 
certain circuit parameters [11, 27] and (2) is the complex 
order parameter (10). The dynamics of x, y and ip are 
given by substitution into (60): 

X — —uy 

y^ux-^-±^{l-x'~y') (62) 
ip = u 

with u = n + Rc{z) - (6 + 1).T. From (56) and (58), 
Re(z) = a;+(-l)"i(a;2 + y2_;^)(^2_^y2)(n-i)/2(,Qgjj^^_ 

(n — 1) ta.n~-^{y/x)]. 

We can now plot the phase portrait for (62) on the 
cyhnder {(a;, y, -0)1 2:, y, -0 G R,x^ + y"^ < 1}. In the 
simple case where a decouples from tp, trajectories can 
be projected down onto the a-disc without intersecting 
themselves or each other. However, in the more typical 
case that a and ■0 are interdependent, we use Poincare 
sections at (mod 27r) = to sort out the structure. In 
these Poincare sections, quasiperiodic trajectories (ide- 
ally) appear as closed curves or island chains, periodic 
trajectories appear as fixed points or period-p points of 
integer period, and chaotic trajectories fill the remaining 
regions of the unit disc. 
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FIG. 2: Poincare sections of a at t/; (mod 27r) — for a resistively-loaded series array of Josephson junctions with b = 1/20, fl = 
3/4 (pt. (A) in Fig. 1). The initial distributions are sinusoidal with wavenumber n, where n is (a) 1, (b) 2, (c) 3, (d) 4, (e) 
5, (f) 6, (g) 7, (h) 8, (i) 16, (j) 32, and (k) oo, i.e. on the Poisson submanifold. In (j) and (k), the complete trajectories are 
plotted instead of the intersections with the plane tp (mod 2tt) = 0. 



First, however, we must choose an appropriate b and 
ri. To do so, we consider their definitions in terms of 
the original circuit parameters: b = R/{NRj) and ft — 
bib/ Ic, where N is the number of junctions, /{, the source 
current, R the load resistance, Ic the critical current of 
each Josephson junction, and Rj the intrinsic Josephson 
junction resistance [11, 27]. Because the resistances must 
be positive in the physical system, we examine only b > Q 
in our simulations. Additionally, Ic represents a positive 
current magnitude, while I^ refiects both a source current 
magnitude and direction. Since the circuit is symmetric 
with respect to reversal of the source circuit (see Fig. 
1 of [27]), the corresponding dynamical system is left 
unchanged by the reflection — > — fi, x — > —x. Hence, 
we also restrict our study to positive values of fi. 

If 6/ri > 1, (62) implies there are fixed points at 
X* — f2/6, y* — ib-y/l — $7^/6^ for arbitrary ip. In nu- 
merical experiments, the negative-y* line of fixed points 
appears to attract distributions, while the positive- j/* line 



repels them. Along the bifurcation curve ^ b, the two 
rows of fixed points merge at a; = 1 , and we find compu- 
tational evidence that a splay state (for which x = y = 0) 
emerges from their union and moves inside the unit disc 
along the a;-axis toward the origin as b is decreased or H. is 
increased. We can see from (62) that any such state must 
lie on the x-axis for all parameter values, as it did in pre- 
vious characterizations of the Poisson submanifold [27]. 

For the submanifolds we examined, chaos only ap- 
peared in the portion of the first quadrant in the b-fl 
plane that did not contain the fixed points, and the chaos 
became more widespread as b/fl — s- 1. This is illus- 
trated schematically in Fig. 1; the gradient of increas- 
ing darkness represents increasingly pervasive chaos. In 
submanifolds where the chaos was not widespread, the 
dynamics on the Poincare sections were reminiscent of a 
Kolmogorov-Arnold-Moser Hamiltonian system with hi- 
erarchies of islands enclosing nested sets of closed orbits. 
Nevertheless, we do not have an explicit Hamiltonian for 
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(a) 



(b) 




(d) 




FIG. 3: Poincare sections of a at i/) (mod 2n) = for a resistively-loaded series array of Josephson junctions with b = 17/10, = 
1 (pt. (B) in Fig. 1). The initial distributions are sinusoidal with wavenumber n, where n is (a) 1, (b) 2, (c) 4, (d) 8, (e) 16, 
(f) 32, (g) 64, and (h) <x, i.e. on the Poisson submanifold. In (g) and (h), the full trajectories are plotted. 



(1) as we do for its averaged counterpart [20]. 

The increase in chaotic behavior is clearly visible in 
Figs. 2 and Fig. 3, which show sequences of Poincare sec- 
tions corresponding to the points (A) and (B) in Fig. 1. 
Point (A) lies at (6, f^) = (1/20,3/4), about 1/2 unit 
from the bifurcation curve — b, while point (B) lies 
at {b,n) = (17/10,1), about 1/3 unit from n ^ b. As 
an example of the pattern of escalating chaos, observe 
that Figs. 3(a),(b),(c) have larger, more dramatically 
overlapping chaotic regions than the corresponding plots 
(a),(b),(d) of Fig. 2. 

Although not shown, the chaotic trajectories that pro- 
duced the scattered points in the Poincare sections are 
phase coherent: they cycle smoothly and unidirection- 
ally around the splay states throughout each period of 
tp. When the splay states are moved toward the edge of 
the unit disc by increasing b or decreasing fi, these tra- 
jectories appear increasingly less prone to return to the 
same neighborhoods in the Poincare sections, resulting 
in the observed amplification of chaotic behavior. 

It is also possible to interpret the association between 
the parameter values and the intensity of the chaos in 
terms of the underlying physical parameters. In terms 
of these parameters, the limit b/il — > 1~ translates to 
Ic/h 1~ or lb J+, which predicts that chaos should 
appear in real series arrays of Josephson junctions if the 
source current is reduced to near the critical current of 
the junctions. 

Even though the Poincare sections in Fig. 2 and Fig. 3 
show differing degrees of chaos, both series of plots de- 
pict a trend of decreasing chaotic behavior with increas- 



ing n. This stems from the dependence of g on the phase 
centroid (z), which in turn arises because the oscillators 
are coupled only through their effect on the first har- 
monic of the phase density. For a coupling of this type, 
a sinusoidal phase density with a very short period and 
rapid oscillations (high n) "looks" nearly identical (in the 
Riemann-Lebesgue sense) to a uniform density. Hence, 
in the limit of large n, we see a decoupling from tjj, just as 
it does on the Poisson submanifold (recall that the Pois- 
son submanifold corresponds to a uniform density in 9, 
as shown in Section IV B). From this perspective, then, 
chaos becomes increasingly dominant as we move "away" 
from the Poisson submanifold, down toward small n. 

Finally, we point out a surprising feature in the 
Poincare sections of (A) that was common in other sim- 
ulations we performed. Starting at n — 5, we see promi- 
nent sets of period- (n — 1) islands which appear for n up 
to 8 in Fig. 2. This ring of islands appears for higher n as 
well and forms an increasingly larger and thinner band 
as n is increased. Inside the dilating band, a set of nested 
orbits resembling the corresponding neutrally stable cy- 
cles of the Poisson submanifold grows, filling the unit disc 
and approaching coincidence with the trajectories on the 
Poisson submanifold. We are currently unclear on why 
exactly (n — 1) islands emerge from the Mobius group ac- 
tion on (58), but pose this as an open question for future 
study. 

Although it is tempting to try to extrapolate our nu- 
merical results to the case of non-identical oscillators, 
Ott and Antonsen [25] have recently demonstrated that 
such systems contain a two-dimensional submanifold (the 
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generalization of the simpler Poisson submanifold stud- 
ied here) that carries all the long-term dynamics of the 
phase centroid (z). Their results hold for the common 
case in which <? is a time-independent angular frequency 
with some distribution of values among the oscillators, 
and / is a function of time, independent of oscillator 
variability. Our numerical experiments, together with 
this new result, indicate that the widespread neutral sta- 
bility in systems of identical, sinusoidally-coupled phase 
oscillators is a consequence of their special symmetries 
and underlying group-theoretic structure. 

VII. APPENDIX 

We show that the 7V!/(iV— 4)! cross ratios of the oscilla- 
tor phases are functionally dependent on the TV — 3 cross 
ratios {A1234, A2345, • • • , A(jv-3)(Ar-2)(Ar-i)Ar}- To do so, 
we use the fact that the 4! cross ratios corresponding 
to the 4! permutations of Zi,Zj,Zk,zi can be written as 
elementary functions of Xijki- 



^ijkl 


= Xjilk = ^klij = ^ 


Xijlk 


= ^/Kjkl 


Mklj 


= 1/(1 — Xijki) 


\kjl 


= 1 — Xijki 


^ilkj 


= Xijki/{1 - Xijki) 


^iljk 


= {Xijki - Xijki 



(63) 



Additionally, we can obtain new cross ratios from existing 
ones by multiplication: 

Xijki Xjmkl = Ximkl (64) 

Using these facts, we need to show that we 
can write Xpqrs for any distinct P,Q,R,S G 
{1,2,..., TV} in terms of elements from 
{^1234, A2345, • • ■ , A(7v-3)(Ar-2)(Af-i)Af}- First, note 
that we can rewrite (64) as a function Fj which takes 
two cross ratios Xijki and Xjkim (with indices in order), 
permutes the indices as necessary to eliminate Zj, 
executes the multiplication and returns the product with 
its indices in order: 

Fj{Xijkl, Xjklm) = Xtklm (65) 



Observe, however, that Fj is just short-hand for a com- 
position of elementary functions from (63): 

Fj{Xijki,Xjkim) = 7 7T TTTT (66) 

J- ~ '^ijkly'^jklm ~ -^jklm 

We can also define the analogous functions and Hi: 

GkiXijkli Xjklm) ~ Xijlm / r.r^\ 

TT(\ \ \-\ ^^^^ 
^ly^ijkli ^jklm) — ^ijkm 

These functions have their own compositions like that of 
Fj in (66). 

Let Xpqrs correspond to the permutation of Xpqrs 
in which the indices are in order. We can write 
XpQRS in terms of Xpqrs using one of the func- 
tions in (63). Thus, the problem reduces to show- 
ing that we can obtain Xpqrs from the elements of 
{Ai234,A2345,---,A(Ar_3)(Ar-2)(JV-i)Ar} by elimination of 
the indices between p, q, r, s using the operations Fj, 
Gk, Hi- 

If there are one or more indices between i and j, we say 

there is a gap between i and j. Now observe that we can 
obtain the first gap between p and q using only Xijki with 
no gaps; we grow this gap iteratively one index at a time 
by the operation: Fk{Xpk{k+i)(k+2),Xk(k+i)(k+2)(k+z)) = 
\{k+i)(k+2)(k+3)- We can then grow the second gap be- 
tween q and r to its full size using only Xiji-i that have 
no gaps between j and k 01 k and I (each of which could 
be made from Xijki with no gaps) using the operation: 
Gk{Xpqk(k+i), Xqk(^k+i)(k+2)) = ^pq{k+i)(k+2)- Finally, we 
can create the third gap between r and s using only 
Xijki with no gaps between k and I (which could be 
made from Xijki with fewer gaps) using the operation: 

Hk{Xpqrk,Xqrk(k+l)) = Xpqr{k+l)- 

Since each Xijki (with i < j < k < I) can be built 
up from Xijki with fewer gaps, the proof is complete: all 
Nl/{N — 4)! cross ratios are dependent on the elements 

of {A1234, A2345, • • • 5 A(Ar-3)(Ar-2)(JV-l)Jv}- 
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